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Abstract 

There  are  many  possible  strategies  for  sampling  the  near-forward  scattering  pattern 
produced  by  a  field  of  large  panicles  and  for  subsequently  solving  the  inverse  problem  in  order  to 
obtain  an  estimate  of  the  particle  size  distribution.  In  a  previous  paper  {Part.  Char.  Vol.  5,  pp.  128- 
133,  1988)  an  optimally  scaled  formulation  of  the  problem  was  derived  based  on  consideration  of 
condition  numbers  of  the  linear  system  obtained  through  numerical  quadrature  of  the  governing 
Fredholm  integral  equation.  Here  we  consider  "scaling"  of  the  problem  to  involve  selection  of  the 
parameters  under  control  of  the  instrument  designer  (e.g.  the  number,  angular  positions,  and 
aperture  geometries  of  the  detectors  and  the  number,  positions,  and  widths  or  weighting  functions 
of  the  discrete  size  classes).  Since  the  many  numericai/analytical  schemes  for  solving  for  the  size 
distribution  given  a  finite  number  of  scattering  measurements  are  fundamentally  very'  similar, 
optimal  scaling  of  the  problem  will  improve  the  performance  of  the  instrument  regardless  of  the 
selected  computational  inversion  algorithm. 

In  this  paper  we  consider  the  analytic  eigenfunction  expansion  method  of  solving  the 
inverse  Fraunhofer  diffraction  problem,  and  in  particular  how  the  scaling  strategy  affects  the 
inversions.  The  eigenfunctions  and  associated  eigenvalues  for  the  diffraction  problem  (assuming 
infinite  support)  are  derived  in  terms  of  two  (variable)  scaling  parameters  which  describe  the 
detector  geometry  and  the  size  class  configurations.  It  is  shown  that  the  rate  of  decrease  of  the 
eigenvalue  spectrum  with  generalized  frequency  is  minimized  for  the  same  scaling  parameters 
(a=2,  b=2)  which  optimized  the  condition  numbers  of  the  linear  system.  (A  slower  decrease  in  the 
eigenvalue  spectrum  indicates  that  the  inversions  will  be  less  susceptible  to  corruption  by  noise,  i.e. 
more  stable).  This  optimal  scaling  study  indicates,  as  was  the  case  for  the  condition  number 
analysis  of  the  discrete  linear  system,  that  the  particle  size  distribution  solution  should  be  on  an  area 
basis  as  n(D)D2  and  the  scattering  measurements  should  provide  i(0)0T 
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In  a  previous  paper  Hirleman  [1]  formulated  the  integral  equation  governing  the  inverse 
Fraunhofer  diffraction  particle  sizing  problem  as  follows: 

oo 

ia(0)  =  /  [  1^X2  /  4*2  ]  [  J,2(a0)  /  (ab-202-a)  ]  nb(a)  da  ( 1 ) 

o 

where:  i(0)  is  the  scattering  intensity  (W/sr)  at  a  small  forward  angle  0;  is  the  irradiance 
(W/m2)  of  the  incident  beam  (assumed  uniform);  a  is  the  nondimensional  particle  size  parameter 
equal  to  the  particle  circumference  divided  by  the  wavelength  X  of  the  incident  radiation  (rtD/X);  D 
is  the  particle  diameter,  n(a)  is  an  unnormalized  particle  frequency  distribution  such  that  n(a)da  is 
the  number  of  optically-sampled  particles  in  the  optical  beam  with  sizes  between  a  and  a  +  da;  J, 
is  a  Bessel  function  of  first  kind  and  first  order,  and  a  and  b  are  constants  where: 

i,(0)  =  i(0)0*  (2) 

n^a)  =  n(a)ab  (3) 

Equation  (1)  assumes  that  Fraunhofer  diffraction  theory  adequately  describes  the  scattered  light 
signature  for  the  angles  0  of  interest  and  that  multiple  scattering  is  negligible.  The  diffraction 
assumption  is  satisfied  generally  when  the  angles  0  are  small,  the  particles  are  large  compared  to 
the  wavelength,  and  the  refractive  index  of  the  particles  relative  to  the  surroundings  is  not  very 
close  to  unity.  Finally,  Eq.  (1)  neglects  coherent  scattering  effects,  i.e.,  assumes  a  very  large 
number  of  randomly  positioned  particles  are  optically  sampled.  Now  the  scaling  of  the  problem, 
which  is  the  under  auspices  of  the  instrument  designer,  involves  selection  of  parameters  a  and  b  and 
the  number  and  orientations  of  detectors  where  measurements  of  i,(0)  will  be  obtained.  The 
measurements  are  used  in  the  inverse  problem,  i.e.  solving  for  nt>(a)  in  Eq.  (1)  using  measured 
i*(0).  Writing  Eq.  (1)  in  a  more  general  integral  equation  form: 

i*(0)  *  ItaeX2/4n?  /  k(a,0)  n^a)  da  (4) 

o 
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where  k(a,8)  is  the  kernel  of  the  integral  equation  and  by  inspection  of  Eq.  (1): 

k(o,8)  =  J,2(a8)  /  (ab-202-a)  (5) 

Now  McWhirter  and  Pike  [21  formulated  a  solution  scheme  to  Fredholm  integral  equations  such  as 
Eqs.  (1)  and  (4)  which  applies  for  cases  where  the  kernel  function  k(a,0)  is  a  function  of  only  the 
product  of  the  two  variables,  i.e.  k  =  k(a9).  With  this  requirement  we  see  from  Eq.  (5)  that: 


2  -  a  =  b  -  2 

(6) 

a  +  b  =  4 

(7) 

and  we  have  lost  one  degree  of  freedom  in  scaling  the  problem  as  a  and  b  are  no  longer 
independent.  Equation  (7)  is  the  same  condition  for  which  the  linear  system  produced  by  numerical 
quadrature  of  Eq.  (1)  takes  the  Toeplitz  form  as  discussed  by  Hirleman  [1].  If  we  define  a  new 
independent  scaling  parameter  5  as: 

8  =  2  -  a  =  b  -  2  (8) 

Then  Eq.  (1)  becomes: 

OO 

*2-5(9)  =  li*^2  /  4tc2  I  kg(a0)  ng+2(a)  da  (9) 

o 

where: 

k§(a0)  =  J12(a0)  /  (a0)5  (10) 

The  eigenfunctions  \|/<o  and  associated  eigenvalues  of  the  kernel,  if  they  exist,  are  defined  by: 

OO 

V®(0)  =  J  k(a0)  da  (11) 

o 

where  to  is  a  continuous  generalized  frequency.  Thus  any  components  of  the  solution  function 
nt»(a)  which  can  be  expressed  in  terms  of  the  eigenfunctions  \|fo»  of  k  will  be  passed  through  the 
integral  operator  intact  but  scaled  by  the  associated  eigenvalue  In  that  case  the  inversion  of  the 
integral  equation  would  be  performed  by  expanding  the  measured  scattering  signature  i2-5(®)  in  a 
series  of  the  eigenfunctions.  The  eigenfunction  components  with  significant  amplitudes  would  then 
be  used  to  synthesize  the  solution  function  after  scaling  by  the  appropriate  eigenvalues.  Now  if  the 
following  integral  of  the  kernel  k  is  bounded: 


/  |k(x)  |x*l/2dx  <  -  (12) 

o 
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then  a  continuum  of  real  eigenfunctions  exist  and  have  been  found  using  Mellin  transforms  by 
McWhirter  and  Pike  [2]  as: 


W(9)  =  Re[(j»o>(e)]  (13) 

ym-(e)=Im[<Me)]  (14) 

where  the  corresponding  real  eigenvalues  are: 

Xra±=  ±|Ar((l/2)  +  io»|  (15) 

and  where  the  Mellin  transform  K(s)  of  k(x)  is  defined  by: 

oo 

K( s)  s  J  Xs'1  k(x)  dx  (16) 

o 

and: 

e-(l/2)-i«>  V  Jt((  1/2)  +  ico) 

W0)  =  — -  (17) 

V(7t|*((l/2)  +  iCO)l) 


It  is  also  sufficient  to  consider  only  (0  >  0  as  pointed  out  by  McWhirter  and  Pike  [2]  and  Viera  and 
Box  [3]. 

Now  we  will  need  to  obtain  the  Mellin  transform  defined  in  Eq.  (16)  of  the  kernel  family  of 
Eq.  (10)  as: 


Kt{ s)  =  /  Xs"1  Ji2(x)  /  x5 
o 

(l/2)5-s  r(&+l-s)  r(l-(5/2)+(s/2)) 
2  H(l+<5/2)-(s/2))  r(2+(5/2)-(s/2)) 


(18) 


We  can  then  calculate  the  eigenvalue  spectrum  from  Eq.  (17)  which  requires: 

(t/2)8-<'/2V'«»r(&Hl/2)-ioj)  r((5/4)-(6/2)+(io>/2)) 

#Ci((l/2)+i<D)  =  -  (19) 

n((3/4H(6/2)-(i<0/2))  T((7/4)+(8/2)-(i(o/2)) 


which,  for  5=2  simplifies  to  Eq.  (3.3)  of  Bertero  and  Pike  [4]. 
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Before  going  further  we  must  also  consider  the  values  of  8  for  which  the  kernel  ks  of  Eq. 
(10)  satisfies  Eq.  (12).  Substituting  k$  for  k  in  Eq.  (12)  results  in  the  condition: 

oo 

!  Jj2(x)  x'(,/2)-8dx  <  <»  (20) 

o 

Convergence  of  this  integral  requires  that  3  >  (-1/2-8)  >  0  or: 

2.5  >  8  >  -0.5  (21) 

For  this  paper  we  consider  integer  values  for  8  of  0,  1,  and  2  which  give  physically  meaningful 
properties  to  the  functions  i2.s  and  ng+2  as  shown  in  Table  I. 


Table  I 

8 

i*»(0) 

Detector 

n6+2(a) 

Form  of  n8+2(a) 

2 

i(0) 

Linear  Array 

a4n(a) 

1 

i(0)0 

Constant  Ar 
Rings 

a^n(a) 

Volume 

distribution 

0 

i(0)02 

Log-scaled 

Rings 

a2n(a) 

Area 

distribution 

Now  an  indication  of  the  quality  of  the  scaling  of  an  inverse  problem  is  the  rate  at  which  the 
eigenvalues  roll  off  with  increasing  frequency  co.  This  is  true  because  the  solution  to  Eq.  (9)  based 
on  eigenfunction  expansions  will  be: 


n**2(a)  =  \  n6aj+  ^(a)  dco  +  J  n8a>-  ^  (a)  d<o  (22) 

o  o 

where: 

ng.aj*  =  (23) 

and: 

W  =  7  i^e)  ^<±(0)  d0  (24) 

o 

It  can  be  seen  from  Eq.  (23)  that  as  the  eigenvalues  A*,*  -*  0,  the  n^co*  and  therefore  the 
solution  grows  without  bound.  The  illconditioned  nature  of  the  problem  would  then  be  manifested 
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as  small  perturbations  in  ij^*,  for  example  due  to  measurement  errors  in  L^O),  would  be 
magnified  greatly  at  frequencies  o>  where  the  are  small.  For  that  reason  it  is  necessary  to 
truncate  the  integrations  in  Eq.  (22)  at  some  finite  value  of  the  generalized  frequency.  In  that 
case  we  approximate  the  solution  as: 


n§+2(a)  =  f  n5co+  dto  +  j  n8(0- ^-(9)  do)  (25) 

o  o 

In  that  case,  components  of  n$f2(<x)  at  frequencies  above  (*)*,„  are  inaccessible  to  the  experiment. 

Now  there  is  a  tradeoff  between  a  desire  to  make  aw*  as  large  as  possible  in  order  to 
minimize  the  truncation  error  in  Eq.  (25)  but  at  the  same  time  keep  the  solution  stable  by  not 
including  frequencies  with  small  eigenvalues.  For  that  reason  the  behavior  of  the  eigenvalue 
spectrum,  shown  in  Fig.  1,  is  very  important.  For  the  possible  values  5  =  0,  1 , 2  the  optimal  value  is 
clearly  8  =  0  which  gives  the  slowest  rolloff  of  The  value  8=0  corresponds  to  a  =  b  =  2,  that 
is  a  solution  on  a  particle  area  basis  using  log-scaled  ring  detectors.  This  result  is  identical  to  that 
obtained  by  Hirleman  [1]  using  a  condition  number  analysis  on  the  linear  system  produced  by 
numerical  quadrature  of  Eq.  (9). 

The  asymptotic  behavior  of  the  eigenvalue  spectrum  is  interesting  as  well,  and  taking  the 
limit  of  Eq.  (19)  as  to  -»  «»  and  using  Eq.  (15)  we  obtain: 

lim  A*,*  =to-«-i  (26) 

<B— »  •• 


Again,  the  value  of  8  =  0  is  the  best  of  the  those  considered,  where  in  the  limit  of  large  to  the 
eigenvalues  are  proportional  to  to-'. 

Some  Comments  on  Instrument  Design 

It  is  clear  from  the  previous  analysis  that  the  parameter  5,  and  in  turn  the  instrument  scaling 
parameters  a  and  b,  control  the  nature  of  the  numerical  inversion  based  on  eigenfunction 
expansions.  Thus  the  stability,  the  information  content,  and  other  measures  of  the  performance  of 
an  inverse  scattering  solution  based  on  analytic  eigenfunctions  is  determined  by  the  selection  of 
parameters  a  and  b.  This  has  also  shown  to  be  the  case  for  inversion  of  the  linear  system  produced 
by  numerical  quadrature  of  Eq.  (9),  and  it  is  expected  to  hold  for  a  singular  function  expansion  [3,5] 
as  well.  It  is  crucial,  then,  to  understand  constraints  on  the  values  of  a  and  b  which  might  be 
implemented  in  practice.  Now  an  underlying  assumption  governing  the  validity  of  the  condition 
number  as  a  measure  of  inversion  stabilityf  1]  and  of  the  scheme  proposed  here  for  truncating  the 
eigenfunction  series  approximation  using  a  noise-based  criterion  is  that  of  white  noise.  If,  for 
example,  the  noise  produced  by  a  detector  were  constant  per  unit  area  of  active  detector  surface, 
then  the  white  noise  assumption  would  only  be  true  for  a  =  0  as  pointed  out  by  Berteno  and  Pike 
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[4,5].  This  forces  the  instrument  (for  a  product  kernel  or  Toeplitz  matrix  form)  to  b  =  4  and,  in 
turn,  8  =  2.  Clearly  it  would  be  better  to  have  the  freedom  to  select  a  detector  geometry  based  on 
optimizing  the  expected  numerical  performance  of  the  inverse  problem  and  force  the  detector  to 
conform  rather  than  vice  versa. 

One  possible  approach  is  to  use  the  optical  system  shown  in  Fig.  2  as  developed  by 
Hirleman  and  Dellenback  [6],  where  the  detector  area  (noise-producing)  does  not  depend  on  the 
collection  aperture  for  the  discrete  detection  angles  0.  The  optical  system  of  Fig.  2  uses  a 
transmission-mode  spatial  light  modulator  to  create  programmable  mosaic  arrays  (of  arbitrary 
shape)  of  detector  openings  which  pass  selected  portions  of  the  scattering  signature  on  through  to 
the  field  lens  and  single  detector.  While  the  parallel,  simultaneous  detection  capabilities  are  lost 
with  a  system  as  in  Fig.  2,  the  potential  for  intelligent,  programmable  detector  arrays  overrides  that 
disadvantage. 
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Nomenclature 


a  instrument  scaling  parameter  (nondimensional)  for  particle  size  distribution  function 

basis 

b  instrument  scaling  parameter  (nondimensional)  for  scattering  intensity  moment 

D  particle  diameter  (m) 

i(0)  scattering  intensity  (W/sr)  at  angle  0 

1^  incident  irradiance  (W/ni2) 

Ji  Bessel  function  of  first  kind  and  first  order 

k(a,0)  general  kernel  function  for  the  integral  equation 

k(a0)  product-form  of  the  kernel  k 

K  Mellin  transform  of  the  kernel  function  k 

a  particle  size  parameter,  rcD/X 

5  instrument  scaling  parameter  (nondimensional)  which  arises  when  k  is  a  product 

kernel. 

0  scattering  angle  measured  from  the  optical  axis 

Xm  eigenvalue  corresponding  to  frequency  to 

eigenfunction  for  the  kernel  of  the  integral  equation 
Q)  generalized  frequency 


log  CO 


Figure  1.  Plot  of  the  eigenvalue  spectrum  X»+  of  the  Fraunhofer  diffraction  integral  equation  for 
the  kernel  described  by  Eq.  (19).  The  spectra  for  three  values  of  the  scaling  parameter  5  are  shown. 
The  asymptotic  behavior  of  the  curves  are  described  by  Eq.  (26). 
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Figure  2.  Schematic  of  a  laser  diffraction  particle  sizing  instrument  where  a  spatial  light  modulator 
operated  as  a  programmable  transmission  mask  has  been  included  to  provide  for  on-line,  adaptive 
configuration  of  the  detector  collection  apertures.  Annular  ring  openings  are  created  in  the  SLM  at 
the  transform  plane  by  setting  the  pixels  to  transmit  or  block  the  incident  polarized  light.  The  rings 
are  created  concentric  with  the  optical  axis  as  measured  in  real-time  by  the  x-y  position  detector 
shown.  The  field  lens  collects  all  light  passing  through  the  SLM  and  passes  it  the  the  detector,  and 
the  system  is  sequenced  through  a  set  of  rings. 
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